function [ varargout ] = divasl_read_asl( varargin )
% Divasl_read_nii Summary of this function goes here
%  If the function is used with path of file(s), then it will read MRI data
%  and then compute the Control-Lable values, as well as the averaged C-L
%  values w.r.t. volumns.
%  If the function is used with data that returned by this same function,
%  then additional parameters are required to specify slices, volumns, TI2
%  ... and another parameter to specify ASL or MRI or ASL_AVERAGE.
%
%  The format of the data is as following:
%  -----------------------------------------------
%   data . mri . data
%              . dim
%        . asl . data
%              . mean
%              . r
%  -----------------------------------------------
%  The method to use it:
%  1. if only file path is at hand, then you can read it by directly using
%         [data_read_from_file] = divasl_read_asl( file_list );
%  2. if the data has been read into the memory, then you can access
%  to every slice (repetitions, TI2) using the following code:
%         [data_sel] = divasl_read_asl(data_read_from_file,
%                                         'Format','asl',...
%                                         'TI2',1,...
%                                         'Slice',8,...
%                                         'Repetition',23);
%  Opions description >>
%          'Format'     -->> which kind of data you want to output
%                        'asl','mri','average'
%          'TI2'        -->> the index of inversion time
%          'Slice'      -->> the index of slices
%          'Repetitions'-->> the index of repetions
%  NOTE: the format of 'data_read_from_file' should be the same as
%  mentioned above.
%
%
%   Dependency: spm_vol, spm_read_vols
%
%   See also: divasl_get_filelist
%

if nargin < 2
    flag = 0;  % flag reseaved for future options, //unused now.//
end

% If reading data from file
if isa(varargin{1},'cell')
    nfiles = length(varargin{1});
    mri = cell(nfiles,1);
    for i = 1:nfiles
        head_vol = spm_vol(varargin{1}{i});
        if iscell(head_vol)
            head_vol = cell2mat(head_vol);
        end
        [data] = spm_read_vols(head_vol);
        mri{i}.data = data;
        [X,Y,Z,R] = size(data);
        mri{i}.dim = struct(...
            'x',X,... % image width
            'y',Y,... % image height
            'z',Z,... % slices
            'r',R); % repetitions
    end
    varargout{1}.mri = mri;
    varargout{1}.asl = asl_label_control(varargout{1},flag);
    return;
end

%--------------------------------------------------------------------------
% If reading data from MRI data read using the same function, additional
% parameters are required to specify fixed slice, TI2, or volumns
if isa(varargin{1},'struct')
    inputs = format_input(varargin{2:end});
    %
    if ~isfield(inputs,'format')
        inputs.format = 'asl';
    end
    switch lower(inputs.format)
        case 'asl'
            data_temp = cell2mat(varargin{1}.asl);
            
            num = length(data_temp);
            dim = size(data_temp(1).data);
            data_out = zeros(dim(1),dim(2),dim(3),dim(4),num);
            for i=1:num
                data_out(:,:,:,:,i) = data_temp(i).data;
            end
        case 'average'
            data_temp = cell2mat(varargin{1}.asl);
            num = length(data_temp);
            dim = size(data_temp(1).mean);
            data_out = zeros(dim(1),dim(2),dim(3),1,num);
            for i = 1:num
                data_out(:,:,:,1,i) = data_temp(i).mean;
            end
        case 'mri'
            data_temp = cell2mat(varargin{1}.mri);
            num = length(data_temp);
            dim = size(data_temp(1).data);
            data_out = zeros(dim(1),dim(2),dim(3),dim(4),num);
            for i = 1:num
                data_out(:,:,:,:,i) = data_temp(i).data;
            end
            
        otherwise
            error('the format you request does not exist.');
    end
    
    dim_sel = cell(ndims(data_out),1);
    
    if isfield(inputs,'slice')
        dim_sel{3} = inputs.slice;
    end

        if isfield(inputs,'ti2')
            if length(dim_sel) <= 3
                disp('Control String TI2 is useless');
            else
                dim_sel{5} = inputs.ti2;
            end
        end
        
        if isfield(inputs,'repetition')
            if length(dim_sel) <= 3
                disp('Control String Repetition is useless');
            else
                dim_sel{4} = inputs.repetition;
            end
        end

    str = '';
    for i = 1:ndims(data_out)
        if isempty(dim_sel{i})
            str = [str,':'];
        else
            str = [str,'[',num2str(dim_sel{i}),']'];
        end
        if i<ndims(data_out)
            str = [str , ','];
        end
    end
    str = ['data_out(',str,')'];
    
    
    varargout{1} = eval(str);
    
    
    return;
end
end

function varargout = format_input(varargin)
num = length(varargin);
if mod(num,2)
    error('input format is not correct');
end
inputs = [];
for i=1:2:num-1
    if ischar(varargin{i})
        inputs = setfield(inputs,lower(varargin{i}),varargin{i+1});
    end
end
varargout{1} = inputs;
end

%--------------------------------------------------------------------------

function varargout = asl_label_control(varargin)
% Basic substraction of control MRI to label MRI
if nargin == 2
    flag = varargin{2};
end
nfiles = length(varargin{1}.mri);
varargout{1} = cell(nfiles,1);
for i = 1:nfiles
    if ~isfield(varargin{1}.mri{i}.dim,'r');
        error(['Error occurs in',mfilename]);
    end
    
    data    = varargin{1}.mri{i}.data;
    volumns = varargin{1}.mri{i}.dim.r;
    vControl = 2; vLabel = 3;
    
    asl_data.r = floor(volumns/2);
    asl_data.data = data(:,:,:,vLabel:2:end,:) - data(:,:,:,vControl:2:end,:);
    asl_data.mean = sum(asl_data.data,4)/asl_data.r;
    
    varargout{1}{i} = asl_data;
end

end